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SUMMARY 

In  1970,  Yefimov  proposed  a  linearized  model  for  predicting  the  effect 
of  water  waves  on  the  adjacent  wind  flow.  We  investigate  the  impact  of 
the  non-linear  terms  neglected  in  his  model.   In  this  investigation,  the 
sea  surface  is  described  by  a  single  harmonic  in  x  and  the  wave  induced 
motion  in  the  atmosphere  is  restricted  to  this  harmonic.   Buoyancy  effects 
are  neglected  and  the  turbulent  fluxes  of  momentum  are  approximated  in  terms 
of  an  eddy  viscosity.  The  equations  are  integrated  numerically  until  a 
steady-state  is  achieved.  These  solutions  show  that  the  nonlinear  terms 
which  were  neglected  by  Yefimov  are  important,  and  that  they  generally 
lead  to  a  smaller  response.  The  numerical  solutions  for  the  mean  wind  are 
compared  with  observations  obtained  by  Davidson  and  Brutsaert.  Although 
the  numerical  solutions  show  some  local  agreement  with  the  observations, 
the  general  behavior  is  at  variance  with  the  observations. 
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Fig.  1.    The  wave  stress,  T  ,  as  a  function  of  z  for  the  linear  model. 

w 

-1  2    -1 

These  curves  are  for  k  =  0.25  m        and  V  =  0.24  m  sec    and  for 

the  following  values  of  U^ :   (1),  0.100  m  sec   ;  (2),  0.175  m  sec   ; 

(3),    0.250  m  sec_1;    (4),    0.275  m  sec"1;    (5),    0.300  m  sec"1; 

(6),    0.350  m  sec"    . 


Fig.  2.    The  wave  stress,  T  ,  as  a  function  of  z  for  the  non-linear  model. 

w 

The  parameters  are  the  same  as  for  Fig.  1. 

Fig.  3.    Isolines  of  u  _,  in  units  of  m  sec   ,  from  Eq.  (38). 

Fig.  4.    Isolines  of  u  _,  in  units  of  m  sec   ,  from  non-linear  numerical 
model  with  K  given  by  (21) . 

Fig.  5.    Same  as  in  Fig.  4  except  using  the  linear  model. 

Fig.  6.    Same  as  in  Fig.  4  except  with  the  K  in  Eq.  (23)  replaced  by  the 
constant  v. 

Fig.  7.    Same  as  in  Fig.  4  except  that  v  varies  with  U.  as  given  by  Eq.  (39) 


1.   Introduction 

In  this  paper  a  numerical  model  is  formulated  to  predict  the  effect  of 
water  waves  on  the  adjacent  air  flow.  Data  described  by  Davidson  (1974)  and 
Brutsaert  (1973)  show  the  effect  of  water  waves  on  the  mean  wind.  Many 
studies  of  the  air  flow  over  waves  have  been  carried  out  for  the  purpose  of 
determining  wave  growth  [Miles  (1957),  Phillips  (1957),  Davis  (1970,  1972), 
Lee  (1972)].  Reynolds  and  Hussain  (1972)  and  Davis  (1970,  1972),  and  Yefimov 
(1970)  emphasized  the  importance  of  including  the  wave-related  turbulent 
Reynolds  stresses;  these  stresses  were  neglected  in  the  linear  model  by 
Miles  (1957).  We  will  follow  the  formulation  of  Yefimov  (1970)  except  that 
we  will  solve  the  initial  value  problem  rather  than  the  steady  state  equations. 
This  is  because  we  wish  to  include  the  nonlinear  interaction  between  the  wave 
induced  motions  and  the  mean  flow.  We  shall  see  that  these  interactions  are 
important . 

This  study  treats  a  single  harmonic,  fully  developed  (not  affected  by 
the  air  motion),  water  wave  field.  The  wave  induced  motion  in  the  atmosphere 
is  also  restricted  to  one  harmonic  in  x  ,  but  the  amplitude  and  phase  of  the 
harmonic  are  free  to  vary.  Actually  as  shown  by  Lee  (1972),  nonlinear 
effects  produce  higher  harmonics.  We  neglect  these  effects  because  we  are 
mainly  concerned  with  the  more  easily  observed  effects  on  the  mean  wind.  The 
equations  for  a  homogeneous  fluid  are  written  in  terms  of  the  Fourier  coeffi- 
cients of  the  wave,  finite  differences  introduced,  and  the  equations  integrated 
until  a  steady  state  is  reached.  These  integrations  are  carried  out  for  a 
number  of  different  cases  and  the  results  are  compared  with  the  observations 
of  Davidson  (1974).   In  order  to  obtain  better  agreement  with  the  observations, 
other  boundary  conditions  and  other  diffusion  coefficients  are  suggested. 
The  nonlinear  effects  on  the  mean  wind  are  investigated  by  comparison  with 
the  linear  solutions. 


Observational  results  used  in  comparisons  were  obtained  from  eddy  correla- 
tion measurements  of  the  momentum  flux  along  with  measurements  of  the  mean 
wind  at  several  levels  and  wave  height  measurements.   The  results  were  inter- 
preted (Davidson,  1974)  with  respect  to  the  dependence  of  the  drag  coefficient 

on  the  non-dimensional  parameter,  c/U,  ,  where  c   corresponds  to  the  phase 

1/2 
speed  at  the  swell  frequency  and   U   is  the  friction  velocity,  U^  =  (f/c) 

This  dependence  of  the  drag  coefficient  on  the  wave  influence  was  extended 

in  order  to  estimate  the  wave  influence  on  the  mean  wind  profile.   The  latter 

is  what  we  examine  in  this  study. 

2.   Basic  equations 

The  equations  for  a  homogeneous  fluid  are 

~  +  v*vV  =   -Vp  +  F,  (1) 

Ot  n~i  n~i  rv 


V   •   V  =  0    ,  (2) 

rv 

where 


.2. 


Fz=S   <-"'"')  +§I   (-«-)•  <*> 


In  this  formulation,  the  turbulent  Reynolds  stresses  are  assumed  to  be  inde- 
pendent of  y  and  p  is  the  non-hydrostatic  pressure  divided  by  the 
density.   The  overbars  represent  local  time  averages,  and  u'   and  w'   are 
the  turbulent  velocity  fluctuations.  The  portion  of  the  velocity  field  which 
is  averaged  over  the  turbulence  time  scale  is  given  by 

V  =  [U(z,t)  +  u(x,z,t)]  i  +  w(x,z,t)  j  ,         (5) 


where  U  is  the  spatially  averaged  (mean)  wind,  and  u  and  w  are  departures 
Thus 

U  =  W  =  0  , 


where   (   )   is  the  average  over  one  wavelength  in  x  . 

We  shall  follow  Yefimov  (1970)  in  representing  the  turbulent  fluxes  by 


dU  ,    ,du  .  d 


w 


-«■«'  -  T-K^  +  v  fe  +  §j)  .  (6) 

2     2 
u,Z  =  w,Z  .  (7) 

Note  that  a  mathematical  treatment  given  later  by  Davis  (1974)  tends  to 
justify  the  eddy  viscosity  form  (6). 

If  (5),  (6)  and  (7)  are  now  substituted  into  (1),  (3)  and  (4),  and  the 

curl  (i  •  Vk)  of  (1)  is  then  taken,  the  following  equation  for  the  vorticity 

.    du   dwx 

(f  =  n—  -  x—)  arises: 

^   dz   dx 

dz  oz  dx  dz   dx 

(Note  in  this  derivation,  dU/dt  has  been  neglected,  and  the  continuity 
equation  (2)  has  also  been  used.) 

Equation  (8)  can  further  be  separated  into  spatially  stationary  and  non- 
stationary  parts  by  averaging  over  one  wavelength,  then  subtracting  the 
averaged  equation  and  neglecting  the  non-stationary  perturbation  products. 
This  yields,  for  the  stationary  part 

dp    d?   _c! /T,du\  /Qv 

udt  +  wdl  =  ~2  (Kd7>  '  (9) 

dz 


and,    for  the  non-stationary  • 

,2„        .2 


^P  ,  t3§  .   d  U     d r  ,du  ,  dw.  ,      r  (3  u   ,  d  w,  ,,_. 

3t  +  U^  +  "^  =  —  [vfe  +  ^)]  "  v[~2-  +  73]    '  (10) 

■  <z    dz  dx  dz   dx 


Since  the  motion  is  non-divergent  (2),  a  stream  function  can  be  defined 
such  that 

u  =  v^  ,   w  =  -  ^*  .  (11) 

oz  ox 


When  these  relations  are  introduced  into  (10),  it  becomes 

ox    oz        dz    dz    dz    ox       dx  dz    dx 


Now  we  introduce  sinusoidal  variations  and  look  for  solutions  which  move 
with  the  speed  c  of  the  water  waves 

i|f  =  A(z,t)  cos  [k(x-ct)]  +  B(z,t)  sin  [k(x-ct)]  .  (13) 

Substitute  (13)  into  (12)  and  separate  the  sine  and  cosine  terms.   Equating 
the  coefficients  of  the  cosine  terms  gives 

jLf^A  _  k2A]    m  k(u.c)(Bk2   .  ^B)+  BkA  +  il[v(^A  +  Ak2)]+  vk2f^A  +  ^    _ 

dz  dz      dz    dz    dz  dz 


(14) 


The  sine  terms  give 


^[^.k2B]=-k(U-c)(Ak2-^)-Ak^  +  ^[V(^+Bk2)]+Vk2[^|-l-Bk2]  .  (15) 

dz  dz     dz    dz    dz  dz 


When  (11)  and  (13)  are  inserted  into  (9)  we  obtain 

dz      oz     dz 

Equations  (14)  -  (16)  clearly  represent  a  coupled,  non-linear  system,  in 

d2U 
that  terms  such  as   B  — t:     occur.   This  essentially  precludes  analytical  solu- 

dz2 

tion  to  the  full  set  of  equations.   Yefimov  proceeded  by  first  considering  the 

functions  involved  to  be  independent  of  time,  thus  reducing  the  system  to  one 
of  ordinary  differential  equations  for  the  steady  state  behavior.   Furthermore, 
in  (14)  -  (15)  he  approximated  the  mean  flow  [our  U(z,t)  term]  by  the  undisturbed  log 
profile,  (U5V/H)  In  (z/z  )  ,  thus  arriving  at  a  linear  system.   After  solving 
this  system,  by  what  he  terms  a  "matrix  distillation"  method,  he  then  solved 
for  the  deviation  from  the  log  profile  by  direct  integration  of  the  steady-state 
counterpart  of  (16).   As  noted  in  the  introduction,  a  major  part  of  our  investi- 
gation concerns  the  effect  of  these  neglected  non-linear  terms.  Also  note,  as 


we  shall  return  to  it  later,  that  (16)  can  be  integrated  twice  to  give  a  simpler 

dU 

)z 


equation  for  ^7  . 


3.   Boundary  Conditions 

Determination  of  the  correct  boundary  conditions  for  (14)  -  (16),  specifically 
at  the  lower  boundary  is  not  as  straightforward  as  the  derivation  of  the 
equations.   Yefimov  (1970)  defined  the  boundary  conditions  on  the  free  surface 
at  the  lower  boundary  to  be 

U(r),t)  +  u(x,T],t)  =  -  akc  cos  k(x-ct)  ,  (17) 

w(x,T],t)  =  akc  sin  k(x-ct)  , 

where  the  height   of  the   free  surface   is  given  by 

r|  =  a  cos  k(x-ct)    .  (18) 


1/2 
Here  a   is  the  wave  amplitude  and  c  =  (g/k)     .   He  then  assumed  that 

the  values  of  U,  u  and  w  on  the  free  surface  can  be  adequately  approxi- 
mated by  their  values  at   z  =  0  ,  if  r)      is  small.   With  u  and  w  replaced 
by  their  values  in  terms  of  (11)  and  (13),  this  reduced  to 


*r;  (0,t)  =  -  akc  ,   A(0,t)  =  ac 


|f  (0,t)  =  0  ,       B(0,t)  =  0  .  (19) 

However,  note  that  (17)  can  be  rewritten 

U(T),t)  +  u(x,n,t)  =  -  ken  , 

hence,  it  seems  more  correct  that  u(x,r|,t)   can  be  approximated  by  its 
value  at  z  =  0  only  if  U(r],t)   can  also  be  (which  is  what  Yefimov  did). 
However,  the  Taylor  expansion  for  U(r),t)   is 

u(n,t)  =  u(o,t)  +  (|2)q  ^  .  (20) 


Since,  for  the  standard  logarithmic  profile,   (t~)^  =  7  »   it  is  very 


,d]L  ~  1 

iZ  o    z 


possible  that  u(x,r],t)  ^u(x,0,t)  .  Note  that  a  similar  difficulty  with 
w(x,0,t)   is  not  likely,  since  U(r),t)  does  not  appear  in  this  condi- 
tion. 

Since  our  main  intent  is  the  investigation  of  the  effect  of  the  non- 
linear terms  in  (14)- (16),  we  shall  use  (19)  as  our  boundary  condition. 
However,  since  our  results  are  at  variance  with  the  observational  studies, 
further  consideration  of  alternative  formulations  to  (19)  seems  warranted. 

Note,  a  consequence  of  Yefimov' s  boundary  conditions  is  that,  near  z  =  0  , 

2 
A(z,t)  =  0(1)   and  B(z,t)  =  0(z  )  . 


Then,  if  one  assumes  that  U(z,t)   is  also  0(1)  ,  as  well  as  second  order 
and  higher  derivatives  of  A  and  B  ,  then  (15)  implies  that 

— ~   =0(1)  ,  near  z  =  0  . 

This  apparently  induced  Yefimov  to  use  — ~  =  0  at  the  lower  boundary  in 

bzZ 
his  numerical  scheme. 

Because  of  the  requirement  that  the  perturbations  die  out  at  the  upper 
boundary   (z  =  H)  ,  both  A(z,t)   and  B(z,t)  ,  as  well  as  their  first 
derivatives  are  set  to  zero  at  z  =  H  . 
4.  Representation  of  the  Mean  Wave  Field 

We  remark  that  in  the  case  of  no  wave  (i.e.,  a  =  0),  the  solution  to 
(14)  -  (16),  plus  the  boundary  conditions,  is  the  log  profile 

A(z,t)  =  B(z,t)  =  0  , 

U* 
U(z,t)  =  —  In  (f-  +   1)  . 

K  Z 

o 
Note,  this  representation  for  the  log  profile  does  not  generate  a  singu- 
larity at  z  =  0  .  Also  note,  from  (16),  this  immediately  implies  a  form 
for  the  diffusion  coefficient 

K  =  x   U,  (z  +  z  )  .  (21) 

o 

Now,  in  the  general  case   (a  ^  0)  ,  there  will  be  a  departure  from 
the  log  profile  which  we  represent 

U* 
U(z,t)  =  —  In  (~  +   1)  +  u(z,t)  ,  (22) 

o 


10 


where 


U^  =  friction  velocity,  z  =  surface  roughness,  h=  von  Karmon's  constant. 

When  (22)  is  substituted  in  equation  (16),  and  the  resulting  relation  integrated, 

we  have 

z 

»--  /     lvdz  <23> 

o 
where  the  wave  induced  stress  is 

T  =  -  uw  =  ■  7  k  T-  A  -  t-  B   .  (24) 

w  2    oz     oz 


We  note  that  Yefimov  (1970)  replaced  the  diffusion  coefficient  K  in  (23)  by 
a  constant  V  ,  which  is  apparently  the  eddy     viscosity  coefficient  used 
in  (6).  Whether  this  was  a  misprint  in  the  translation  is  not  clear,  and 
we  investigated   the  solutions  for  both  K  as  given  by  (21)  and  with  K  =  v  . 
5.  Numerical  Solution  Procedure 

As  noted  above,  Yefimov  solved  the  linearized  steady-state  equivalents  of 
(14)  -  (15)  by  a  procedure  described  only  by  the  term  "matrix  distillation." 
In  contrast,  we  chose  to  solve  the  time-dependent  non-linear  system  (14),  (15), 
(22)  -  (24)  by  difference  equations,  and  time  iterate  until  a  numerical  steady 
state  evolved.   To  implement  this,  (14)  and  (15)  were  first  reduced  to  difference 
equation  form,  using  central  finite  differences  on  all  spatial  operations  and  all 
time  operations  except  the  diffusion  terms  (i.e.,  those  with  coefficient  v  )• 
The  latter  were  evaluated  at  the  previous  time  step  in  order  to  eliminate  a 
diffusive  instability.   The  left  hand  sides  of  these  equations  were  then  solved 
for  the  true  tendencies  with  the  Gauss  elimination  technique  described 


11 


in  Richtmyer  (1957).   Leapfrog  time  differencing  was  used  except  that  the 
integration  was  restarted  every  50  time  steps  with  the  finite  difference 
scheme  developed  by  Matsuno  (1966).   This  procedure  eliminates  solution 
separation. 

The  values  of  A(z,t)  and  B(z,t)  thus  computed  are  used  to  update  the 
mean  wind,  as  given  by  (21)  -  (24)  ,  and  this  updated  mean  wind  is  then  used 
to  compute  A(z,t)  and  B(z,t)  at  the  next  time  step.  The  boundary  condi- 
tions take  the  following  finite  difference  form 

A(0,t)  =  ac  ,  (25) 

B(0,t)  =  0  ,  (26) 


A(Az,t)  -  A(0,t)  m  m   akc 

Az  ' 

B(Az,t)  -  B(0,t)  m   Q 

Az 


(27) 


(28) 

A(H,t)  =  B(H,t)  =  0  ,  (29) 

A(H.t)-  A(H-Az.t)  m   B(H,t)-  B(H-At)  = 

Az  Az  '  v 

where     Az      is   the  grid   increment. 

The   following   initial  conditions  were  used  to  start   the  solution 

A(z,0)   =  ac   sinh    [MIL-   z)]/sinh  kH      ,  (31) 

B(z,0)   =   0    , 

U(z,0)   =  —  In   (-5-  +  1)    ,  (32) 

K  Z 

o 

where  H.  =  H-Az  .   Note,  the  fields  (31)  represent  potential  flow  over  the 
waves . 

From  these  initial  conditions,  as  noted  above,  the  time  integration  of 
the  system  of  equations  is  accomplished  by  first  predicting  A  and  B   from 
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(14)  and  (15)  and  then  by  obtaining  U  from  (21)  -  (24).   The  integration  is 

continued  until  the  quantity 

f  "    1 

AZ  2   1/2 

ST=    [    Z        [T(zi,t)   -   T(z.,t  -   At)p]i/Z    ,  (33) 

i=l 

satisfies   the  condition 

ST  <  10"5m2sec_1    .  (34) 

Condition  (34)  is  taken  to  describe  steady  state. 
6.   Observational  studies 

Our  numerical  results  will  be  compared  to  the  results  of  several  observa- 
tional studies  which  have  considered  this  problem  and  have  looked  for  a 
dependence  on  c/U^.  .   Brutsaert  (1973)  has  suggested  that  the  wind  shear 
should  have  the  following  dependence 

i=|[l  +  P(^-a)2]  ,  (35) 

for  c/U,  >  a  •   From  the  data  he  found  that 

(3  =  .006  and  a  =  29  . 

Davidson  (1974)  examined  the  wind  at  10  meters  with  the  BOMEX  data  and 
obtained  the  following  relation 

uio  =.4T  [ln(r>  +  6"57  I  +  •"&  '  26'3)1  •       <36) 

O  * 


These  studies  show  increased  wind  due  to  the  wave  for  c/U^  >  a  where  a  is 
between  25  and  30.   We  will  see  later  that  the  assumption  that  u   is  dependent 
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on  only  c/U^  is  probably  too  restrictive.   In  particular  Brutsaert  (1973) 
assumed  that  the  wave  stress  is  constant  in  height  which  is  neither  reasonable 
nor  consistent  with  observations. 
7.   Numerical  results 

The  numerical  model,  as  described  above,  was  integrated  for  a  number  of 
the  waves  considered  by  Yefimov  (1970)  and,  as  a  test,  the  system  was  made 
linear  by  setting  u  equal  to  zero.  When  the  dimensions  were  restored  to 
Yefimov 's  solution,  the  comparison  was  quite  good.   It  should  be  noted  here 
that  Yefimov  non-dimensionalized  his  variables  by  scaling,  and  the  results  he 
presented  corresponded  to  unsealed  values  of  U^  in  the  range  0.01  msec   < 
U^  <  0.05  msec   ,   or  approximately  an  order  of  magnitude  smaller  than  the 
values  for  which  Davidson  and  Brutsaert  obtained  experimental  data.  The  results 
which  will  be  shown  use  the  following  values  for  the  various  parameters,  unless 
otherwise  noted: 

AZ  =  0.25m 

At  o  0.0125sec 

H,  =  20m 

k  =  0.35 

z  =  3  x  10~4m  (37) 

o 

a  =  0.1/k 

The  value  for  H-   was  selected  after  numerical  experiments  showed  larger 
values  did  not  significantly  alter  the  wind  field  at  the  values  of  c  and  U^ 
of  interest. 

Fig.  1  shows  a  typical  stress  field  prediction  for  the  linear  (Yefimov) 
model  for  different  values  of  U^,  at  wave  number  k  =  0.25m    and 
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2   -1 
V  =  0.24  m  sec    .   The  values  for  k  and  U   are  typical  for  the  BOMEX 

data  used  by  Davidson,  and  the  value  for  v  is  that  used  by  Yefimov.   All 
of  the  curves  have  the  same  general  pattern,  with  initially  growing  ampli- 
tude with  increasing  UA,  up  to  some  maximum  amplitude,  followed  by  decreas- 
ing amplitude  for  further  increases  in  U^.   Yefimov  referred  to  this 
as  a  resonance-like  behavior. 

Fig.  2  presents  the  stress  field  predicted  for  the  same  wave  number  and 
values  of  U^  using  the  full  non-linear  model.   Observe  the  qualitative 
behaviors  are  similar;  however,  two  important  quantitative  differences  emerge. 
First,  the  magnitude  of  the  maximum  stress  amplitude  for  the  non-linear  model 
is  only  about  35^  of  that  predicted  by  the  linear  model.  We  found  this 
pattern  consistently  through  all  ranges  of  wave  numbers  considered,  with 
the  maximum  stress  amplitude  predicted  by  the  non-linear  model  usually  only 
one-quarter  to  one-half  that  of  the  linear  model.   Secondly,  note  the  apparent 
"resonance"  is  noticeably  less  pronounced  in  the  non- linear  case.  This 
pattern  also  appeared  consistently  over  the  range  of  wave  numbers.   From 
these  observations,  we  have  concluded  that  the  presence  of  the  non-linear 
terms  in  (14)  -  (16)  acts  strongly  to  damp  the  wave  induced  stresses,  and 
that  to  approximate  (14)  -  (16)  linearly,  as  Yefimov  did,  is  not  valid. 
Some  further  calculations  were  carried  out  with  ak  =  .05  rather  than  the 
value  ak  =  0.1  that  was  used  by  Yefimov.  With  this  relation  for  a  the 
comparison  between  the  linear  and  nonlinear  solutions  was  closer  except 
for  the  longer  waves  near  k  =  .05  m   .  As  a  final  point,  note  that  both 
models  predict  the  main  contributions  to  the  stress  field  occur  in  the  lowest 
5-8  meters.   This  is  especially  important,  since  both  experimental  studies 
referred  to  above  have  limited  or  non-existent  data  in  this  region. 
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In  comparing  these  figures  to  the  wind  shear  hypothesized  by  Brutsaert, 
(35),  it  is  obvious  that  neither  the  linear  nor  the  non-linear  models  agree 
with  his  assumption  of  constant  stress  at  all  heights.  As  noted  above,  we 
feel  this  assumption  was  too  restrictive,  though  until  additional  data  is 
obtained  in  the  region  below  8  meters  it  cannot  be  fully  discounted. 

To  compare  the  results  of  the  numerical  models  with  Davidson's,  we  ex- 
amined the  deviation  from  the  logarithmic  profile  predicted  by  Davidson: 

"u>  =  uio  -  o74i  ln  (r>  ■  §31  (c  "  26'3  D*>  '  (38) 

o 

c  >  27.5  UA 

for  the  neutral  case,  with  our  predicted  steady  state  value  of  u(10,t)  , 
over  a  range  of  representative  values  of  c  and  U*  .   Fig.  3  shows  the  level 
curves  of  u.ft  predicted  by  Davidson.  Note  that,  even  allowing  for  uncer- 
tainties in  the  determination  of  numerical  coefficients,  Davidson  predicts 
essentially  linear  isolines.   Fig.  4  shows  the  isolines  determined  by  our 
non-linear  model,  using  (21)  to  determine  the  diffusion  coefficient,  K  . 
Although  there  are  obviously  regions,  especially  near  the  line  c  =  27.5  UA  , 
where  locally  linear  isolines  occur,  it  is  clear  that  our  model  will  not 
produce  this  behavior  globally.  More  seriously,  although  generally  our 
model  predicts  increasing  magnitude  for  the  deviation  as   (c/UA)   increases, 
which  is  consistent  with  Davidson's  observations,  the  algebraic  sign  is  not 
consistent  for  observations  show  a  positive  deviation  growing  in  magnitude, 
while  our  model  predicts  a  negative  deviation  growing  in  magnitude. 

The  observed  sign  of  u  at  10  meters  could  be  obtained  by  reversing 
the  curves  in  Fig.  2  or  by  distorting  them  in  such  a  way  that  the  negative 
area  would  be  larger  than  the  positive  area.  Davidson  (1974)  did  observe 
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negative  stresses  at  the  3  meter  level  and  in  an  earlier  study  with  Lake 
Michigan  data  Davidson  and  Frank  (1973)  found  positive  stresses  in  the 
lower  layers.   It  is  not  clear  at  present  how  our  model  could  be  restruc- 
tured to  reverse  this  algebraic  sign,  although  we  suspect  that  the 
boundary  conditions  are  the  major  influence.   There  are  two  specific  areas 
of  the  boundary  conditions  that  we  feel  need  more  study.   The  first  of 
these  is  the  adequacy  of  replacing  values  for  the  horizontal  components 
on  the  free  surface  by  their  values  on  the  mean  surface,  as  was  discussed 
in  section  3.   A  second  point  is  the  validity  of  the  parameterization  of 
the  turbulent  stresses  used  by  Yefimov  (see  Eqs.  (6)  and  (7)).   The  be- 
havior of  these  stresses  near  the  water  surface  is  particularly  crucial. 

Figs.  5,  6,  and  7  show  the  isolines  predicted  by  other  variants  on 
the  model.   Those  in  Fig.  5  are  derived  from  the  linear  model.   Those  in 
Fig.  6  are  obtained  from  our  non-linear  model, by  replacing  the  variable 
diffusion  coefficient,  K  in  (23),  by  the  constant  v  ,  as  Yefimov's  paper 
indicated.   The  curves  in  Fig.  7  are  derived  using  a  V  which  varied 
with  U.   according  to 

v  =   1.2  sec"1  U^  .  (39) 

This  was  suggested  by  the  development  of  Davis  (1974).   Except  for  the 
linear  model's  unacceptably  large  predictions  for  the  deviation  (due  to  the 
"damping"  we  hypothesized  in  the  non-linear  model) ,  all  of  these  share  the 
same  qualitative  features  as  we  discussed  for  Fig.  4.   This  essentially 
leads  us  to  believe,  as  mentioned  above,  that  any  changes  in  algebraic  sign 
of  the  deviations  must  be  through  changes  in  the  boundary  conditions,  or 
perhaps  through  the  formulation  of  the  wave  induced  turbulent  diffusion. 
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8.   Conclusions 

The  purpose  of  this  paper  was  to  investigate  the  linear  model  proposed 
by  Yefimov  from  two  aspects.   First,  to  determine  the  effect  of  the  non- 
linear terms  neglected  by  him,  and  secondly  to  determine  the  agreement 
between  predictions  by  the  linear  and  non-linear  models  and  observational 
studies  by  Davidson  (1974)  and  Brutsaert  (1973),  which  appeared  after 
Yefimov's  work. 

The  main  conclusions  of  this  study  are  that  the  non-linear  terms  neg- 
lected by  Yefimov  are  in  fact  quite  important  when  ak  =  0.1  ,  and  impart 
a  strong  damping  to  the  induced  fields,  that  neither  a  linear  nor  a  non- 
linear model  will  produce  a  height- independent  stress  as  predicted  by 
Brutsaert,  and  that  although  the  non- linear  model  produces  some  agreement 
with  Davidson's  observations,  there  are  also  important  disagreements,  the 
most  important  being  in  algebraic  sign.  However,  since  the  model  predicts 
the  major  contribution  to  the  induced  stress  to  occur  in  the  region  near 
the  mean  surface,  and  since  there  are  significant  questions  about  the 
proper  type  of  the  mathematical  boundary  conditions,  resolution  of  the 
differences  and  determination  of  how  the  model  should  be  modified  to  pro- 
duce results  more  consistent  with  observation  demands  more  complete  obser- 
vational data  on  the  induced  stresses  near  the  surface. 
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